【編集履歴】

  • 2023/10/18:執筆開始
  • 2023/10/19:「アルキメデスの螺旋の可視化」を追加

アルキメデスの螺旋の可視化

 sin関数(サイン関数・sine function)とcos関数(コサイン関数・cosine function)を用いて定義されるアルキメデスの螺旋(アルキメデスの渦巻・Archimedes’ spiral)をグラフで確認します。

 利用するパッケージを読み込みます。

# 利用パッケージ
library(tidyverse)
library(gganimate)

 この記事では、基本的に パッケージ名::関数名() の記法を使うので、パッケージの読み込みは不要です。ただし、作図コードについてはパッケージ名を省略するので、 ggplot2 を読み込む必要があります。
 また、ネイティブパイプ演算子 |> を使います。magrittr パッケージのパイプ演算子 %>% に置き換えられますが、その場合は magrittr を読み込む必要があります。

定義式の確認

 まずは、アルキメデスの螺旋の定義式を確認します。
 sin関数については「sin関数の可視化」、cos関数については「cos関数の可視化」、度数法と弧度法の角度については「円周の作図」を参照してください。

 アルキメデスの螺旋は、sin関数とcos関数を用いて、次の式で定義されます。

\[ \begin{align*} &\quad r = a \theta \\ &\left\{ \begin{aligned} x &= r \cos \theta \\ y &= r \sin \theta \end{aligned} \right. \end{align*} \]

 変数 \(\theta\) は弧度法における角度(ラジアン)です。度数法における角度を \(\theta^{\circ}\) とすると、\(\theta = \frac{2 \pi}{360^{\circ}} \theta^{\circ}\) の関係です。\(\pi\) は円周率です。
 螺旋の形状は定数(パラメータ)によって決まります。実数 \(a\) は螺旋の間隔や全体のサイズに影響します。
 \(r\) は中心からの距離で、\(\theta\) に応じて変化します。

アルキメデスの螺旋の作図

 次は、アルキメデスの螺旋のグラフを作成して、変数(ラジアン)と座標(曲線上の点)の関係を確認します。

螺旋の作図

 アルキメデスの螺旋のグラフを作成します。

 パラメータを指定して、螺旋の描画用のデータフレームを作成します。

# パラメータを指定
a <- 1

# 周回数を指定
lap_num <- 4

# 螺旋用のラジアン(弧度法の角度)を作成
t_vec <- seq(from = 0, to = lap_num*2*pi, length.out = 1000) # 正の範囲
#t_vec <- seq(from = -lap_num*2*pi, to = 0, length.out = 1000) # 負の範囲
#t_vec <- seq(from = -lap_num*pi, to = lap_num*pi, length.out = 1000) # 正と負の範囲

# 螺旋の座標を作成
spiral_df <- tibble::tibble(
  t = t_vec, # ラジアン
  r = a * t, # ノルム
  x = r * cos(t), # x座標
  y = r * sin(t), # y座標
  sgn_t_flag = t >= 0 # 角度の正負
)
spiral_df
## # A tibble: 1,000 × 5
##         t      r      x        y sgn_t_flag
##     <dbl>  <dbl>  <dbl>    <dbl> <lgl>     
##  1 0      0      0      0        TRUE      
##  2 0.0252 0.0252 0.0251 0.000633 TRUE      
##  3 0.0503 0.0503 0.0503 0.00253  TRUE      
##  4 0.0755 0.0755 0.0753 0.00569  TRUE      
##  5 0.101  0.101  0.100  0.0101   TRUE      
##  6 0.126  0.126  0.125  0.0158   TRUE      
##  7 0.151  0.151  0.149  0.0227   TRUE      
##  8 0.176  0.176  0.173  0.0309   TRUE      
##  9 0.201  0.201  0.197  0.0402   TRUE      
## 10 0.226  0.226  0.221  0.0508   TRUE      
## # ℹ 990 more rows

 周回数 lap_num を指定して、螺旋の座標計算用のラジアン \(\theta\) を作成します。\(2 \pi\) の範囲で1周します。
 実数 \(a\) を指定して、原点からの距離(変数ごとの係数) \(r = a \theta\) を計算し、x軸の値 \(x = r \cos \theta\) とy軸の値 \(y = r \sin \theta\) を計算します。

 グラフサイズや軸線用の値を設定します。

# グラフサイズを設定
axis_size <- c(spiral_df[["x"]], spiral_df[["y"]]) |> 
  abs() |> 
  max() |> 
  ceiling()

# 目盛間隔を設定
step_val <- 15

# 軸線数を設定
if(axis_size%%step_val == 0) {
  
  # グラフサイズと最大目盛が一致する場合
  circle_num <- axis_size %/% step_val
} else {
  
  # グラフサイズと最大目盛が一致しない場合
  circle_num <- axis_size %/% step_val + 1
  
  # グラフサイズを拡大
  axis_size <- max(axis_size, circle_num*step_val)
}
axis_size; circle_num
## [1] 30
## [1] 2

 ノルム軸線の目盛間隔(円ごとの半径の間隔)を step_val、軸線数(円の数)を circle_num とします。作図時に自動で設定されるx軸・y軸の目盛間隔と同じ値にすると見栄えがよくなります。
 螺旋のx軸・y軸の最大値とノルム軸の最大値からグラフサイズの半分を axis_size とします。

 ノルム軸線(円)の描画用のデータフレームを作成します。

# ノルム軸線の座標を作成
coord_circle_df <- tidyr::expand_grid(
  r = 1:circle_num * step_val, 
  t = seq(from = 0, to = 2*pi, length.out = 361), 
) |> # 半径ごとにラジアンを複製
  dplyr::mutate(
    x = r * cos(t), 
    y = r * sin(t)
  )
coord_circle_df
## # A tibble: 722 × 4
##        r      t     x     y
##    <dbl>  <dbl> <dbl> <dbl>
##  1    15 0       15   0    
##  2    15 0.0175  15.0 0.262
##  3    15 0.0349  15.0 0.523
##  4    15 0.0524  15.0 0.785
##  5    15 0.0698  15.0 1.05 
##  6    15 0.0873  14.9 1.31 
##  7    15 0.105   14.9 1.57 
##  8    15 0.122   14.9 1.83 
##  9    15 0.140   14.9 2.09 
## 10    15 0.157   14.8 2.35 
## # ℹ 712 more rows

 半径が step_val 間隔の circle_num 個の円周の座標を計算します。外側の円周の半径は axis_size になります。

 角度軸線(斜線)の描画用のデータフレームを作成します。

# 半円における目盛数(分母の値)を指定
denom <- 6

# 角度軸線の座標を作成
coord_oblique_df <- tibble::tibble(
  i = seq(from = 0, to = 2*denom-1, by = 1), # 目盛位置番号(分子の値)
  t = i / denom * pi, 
  r = axis_size, # ノルム軸線の最大値
  x = r * cos(t), 
  y = r * sin(t), 
  t_label = paste0("frac(", i, ", ", denom, ")~pi"), # 角度ラベル
  h = 1 - (x/r * 0.5 + 0.5), 
  v = 1 - (y/r * 0.5 + 0.5)
)
coord_oblique_df
## # A tibble: 12 × 8
##        i     t     r         x         y t_label             h      v
##    <dbl> <dbl> <dbl>     <dbl>     <dbl> <chr>           <dbl>  <dbl>
##  1     0 0        30  3   e+ 1  0        frac(0, 6)~pi  0      0.5   
##  2     1 0.524    30  2.60e+ 1  1.5 e+ 1 frac(1, 6)~pi  0.0670 0.25  
##  3     2 1.05     30  1.5 e+ 1  2.60e+ 1 frac(2, 6)~pi  0.25   0.0670
##  4     3 1.57     30  1.84e-15  3   e+ 1 frac(3, 6)~pi  0.5    0     
##  5     4 2.09     30 -1.5 e+ 1  2.60e+ 1 frac(4, 6)~pi  0.75   0.0670
##  6     5 2.62     30 -2.60e+ 1  1.5 e+ 1 frac(5, 6)~pi  0.933  0.25  
##  7     6 3.14     30 -3   e+ 1  3.67e-15 frac(6, 6)~pi  1      0.5   
##  8     7 3.67     30 -2.60e+ 1 -1.5 e+ 1 frac(7, 6)~pi  0.933  0.75  
##  9     8 4.19     30 -1.50e+ 1 -2.60e+ 1 frac(8, 6)~pi  0.75   0.933 
## 10     9 4.71     30 -5.51e-15 -3   e+ 1 frac(9, 6)~pi  0.5    1     
## 11    10 5.24     30  1.5 e+ 1 -2.60e+ 1 frac(10, 6)~pi 0.25   0.933 
## 12    11 5.76     30  2.60e+ 1 -1.50e+ 1 frac(11, 6)~pi 0.0670 0.75

 ノルム軸線の外側の円周上を等間隔の 2 * denom 個に分割した点の座標を計算します。denom\(n\) とすると、目盛間隔は \(\frac{1}{n} \pi\) になります。また \(i, n\) を整数として、\(\frac{i}{n} \pi\)\(\frac{i + 2 n}{n} \pi = \frac{i}{n} \pi + 2 \pi\) の目盛位置が一致します。

 螺旋のグラフを作成します。

# ラベル用の文字列を作成
fnc_label <- paste0(
  "list(", 
  "r == a * theta, ", 
  "a == ", a, 
  ")"
)

# 螺旋を作図
add_size <- 2
ggplot() + 
  geom_path(data = coord_circle_df, 
            mapping = aes(x = x, y = y, group = r), 
            color = "white") + # ノルム軸線
  geom_segment(data = coord_oblique_df, 
               mapping = aes(x = 0, y = 0, xend = x, yend = y, group = i), 
               color = "white") + # 角度軸線
  geom_text(data = coord_oblique_df, 
            mapping = aes(x = x, y = y, label = t_label, hjust = h, vjust = v), 
            parse = TRUE) + # 角度目盛ラベル
  geom_path(data = spiral_df, 
            mapping = aes(x = x, y = y), 
            linewidth = 1) + # 螺旋
  coord_fixed(ratio = 1, 
              xlim = c(-axis_size-add_size, axis_size+add_size), 
              ylim = c(-axis_size-add_size, axis_size+add_size)) + 
  labs(title = "Archimedes' spiral", 
       subtitle = parse(text = fnc_label), 
       x = expression(x == r ~ cos~theta), 
       y = expression(y == r ~ sin~theta))

 (x軸方向に値が増減するため geom_line() ではなく、) geom_path() で螺旋を描画します。

 \(2 \pi\) の間隔で1周するので、線の間隔は \(2 a \pi\) で一定になります。

 以上でアルキメデスの螺旋を描画できました。次からは、変数やパラメータによる螺旋への影響を確認していきます。

角度と座標の関係

 アルキメデスの螺旋上を移動する点のアニメーションを作成します。

 フレーム数を指定して、螺旋上の点の描画用のデータフレームを作成します。

# フレーム数を指定
frame_num <- 300

# 螺旋上の点用のラジアンを作成
t_vals <- seq(from = min(t_vec), to = max(t_vec), length.out = frame_num+1)[1:frame_num]

# 螺旋上の点の座標を作成
angle_point_df <- tibble::tibble(
  frame_i = 1:frame_num, # フレーム番号
  t = t_vals, 
  r = a * t, 
  x = r * cos(t), 
  y = r * sin(t), 
  var_label = paste0(
    "list(", 
    "a == ", a, ", ", 
    "theta == ", round(t/pi, digits = 2), " * pi, ", 
    "r == ", round(r, digits = 2), 
    ")"
  ) # 変数ラベル
)
angle_point_df
## # A tibble: 300 × 6
##    frame_i      t      r      x       y var_label                               
##      <int>  <dbl>  <dbl>  <dbl>   <dbl> <chr>                                   
##  1       1 0      0      0      0       list(a == 1, theta == 0 * pi, r == 0)   
##  2       2 0.0838 0.0838 0.0835 0.00701 list(a == 1, theta == 0.03 * pi, r == 0…
##  3       3 0.168  0.168  0.165  0.0279  list(a == 1, theta == 0.05 * pi, r == 0…
##  4       4 0.251  0.251  0.243  0.0625  list(a == 1, theta == 0.08 * pi, r == 0…
##  5       5 0.335  0.335  0.316  0.110   list(a == 1, theta == 0.11 * pi, r == 0…
##  6       6 0.419  0.419  0.383  0.170   list(a == 1, theta == 0.13 * pi, r == 0…
##  7       7 0.503  0.503  0.440  0.242   list(a == 1, theta == 0.16 * pi, r == 0…
##  8       8 0.586  0.586  0.488  0.325   list(a == 1, theta == 0.19 * pi, r == 0…
##  9       9 0.670  0.670  0.525  0.416   list(a == 1, theta == 0.21 * pi, r == 0…
## 10      10 0.754  0.754  0.550  0.516   list(a == 1, theta == 0.24 * pi, r == 0…
## # ℹ 290 more rows

 フレーム数 frame_num を指定して、螺旋の座標計算用のラジアン t_vec の範囲を等間隔に frame_num 個に分割して、螺旋上の点の座標計算用のラジアン t_vals を作成します。

 角度を示す線分の描画用のデータフレームを作成します。

# 反転フラグを設定
neg_flag <- FALSE

# 角度線の座標を作成
angle_oblique_df <- dplyr::bind_rows(
  # 角度用の斜線
  tibble::tibble(
    frame_i = 1:frame_num, 
    t = t_vals, 
    r = a * t, 
    x = dplyr::if_else(r >= 0, true = axis_size*cos(t), false = -axis_size*cos(t)), 
    y = dplyr::if_else(r >= 0, true = axis_size*sin(t), false = -axis_size*sin(t)), 
    type = "origin"
  ), 
  # 反転した角度用の斜線
  tibble::tibble(
    frame_i = 1:frame_num, 
    t = t_vals, 
    r = a * t, 
    x = dplyr::case_when(
      neg_flag == FALSE ~ NA, 
      r <  0 ~ axis_size * cos(t), 
      r >= 0 ~ -axis_size * cos(t)
    ), 
    y = dplyr::case_when(
      neg_flag == FALSE ~ NA, 
      r <  0 ~ axis_size * sin(t), 
      r >= 0 ~ -axis_size * sin(t)
    ), 
    type = "negation"
  )
)
angle_oblique_df
## # A tibble: 600 × 6
##    frame_i      t      r     x     y type  
##      <int>  <dbl>  <dbl> <dbl> <dbl> <chr> 
##  1       1 0      0       30    0    origin
##  2       2 0.0838 0.0838  29.9  2.51 origin
##  3       3 0.168  0.168   29.6  5.00 origin
##  4       4 0.251  0.251   29.1  7.46 origin
##  5       5 0.335  0.335   28.3  9.87 origin
##  6       6 0.419  0.419   27.4 12.2  origin
##  7       7 0.503  0.503   26.3 14.5  origin
##  8       8 0.586  0.586   25.0 16.6  origin
##  9       9 0.670  0.670   23.5 18.6  origin
## 10      10 0.754  0.754   21.9 20.5  origin
## # ℹ 590 more rows

 \(\theta\) に応じて、半径が axis_size の円周上の点の座標を計算します。\(r\) が負の値の場合は、螺旋が反転するので、半径を -axis_size とします。
 \(r\) が正負の値を含む場合は、反転した点の座標も格納します。(全ての場合を共通のコードで処理するために、)反転した座標を含めるかどうかを neg_flag で設定します。

 螺旋上の点の描画用のデータフレームを作成します。

# 角度線上の点の座標を作成
r_min <- a * min(t_vec)
r_max <- a * max(t_vec)
point_df <- tidyr::expand_grid(
  frame_i = 1:frame_num, 
  lap_i   = (-lap_num):(lap_num-1) # 点番号
) |> # フレームごとに点を複製
  dplyr::mutate(
    tmp_t = t_vals[frame_i] %% (2*pi), # 単位円相当のラジアン
    t = lap_i * 2*pi + tmp_t, 
    r = a * t, 
    x = r * cos(t), 
    y = r * sin(t),   
    sgn_t_flag = t >= 0 # 角度の正負
  ) |> 
  dplyr::filter(r >= r_min, r <= r_max) # 螺旋上の点を抽出
point_df
## # A tibble: 1,200 × 8
##    frame_i lap_i  tmp_t       t       r       x         y sgn_t_flag
##      <int> <int>  <dbl>   <dbl>   <dbl>   <dbl>     <dbl> <lgl>     
##  1       1     0 0       0       0       0       0        TRUE      
##  2       1     1 0       6.28    6.28    6.28   -1.54e-15 TRUE      
##  3       1     2 0      12.6    12.6    12.6    -6.16e-15 TRUE      
##  4       1     3 0      18.8    18.8    18.8    -1.38e-14 TRUE      
##  5       2     0 0.0838  0.0838  0.0838  0.0835  7.01e- 3 TRUE      
##  6       2     1 0.0838  6.37    6.37    6.34    5.33e- 1 TRUE      
##  7       2     2 0.0838 12.7    12.7    12.6     1.06e+ 0 TRUE      
##  8       2     3 0.0838 18.9    18.9    18.9     1.58e+ 0 TRUE      
##  9       3     0 0.168   0.168   0.168   0.165   2.79e- 2 TRUE      
## 10       3     1 0.168   6.45    6.45    6.36    1.08e+ 0 TRUE      
## # ℹ 1,190 more rows

 \(\theta\) の前後 \(2 \pi\) 間隔 \(\theta + 2 i \pi\) の点の座標を計算します。\(i\) は整数で、lap_i 列とします。
 lap_i 列は、\(r\) が正の値のみの場合は 0 から lap_num - 1、負の値のみの場合は -lap_num から -1、正負の値を含む場合は ceiling(-0.5 * lap_num) から floor(0.5 * lap_num - 1) の整数です。(全ての場合を共通のコードで処理するために、)-lap_num から lap_num - 1 の整数を用いて座標を計算し、r の最小値から最大値の点(行)を取り出します。

 螺旋上の点のアニメーションを作成します。

# 螺旋上の点のアニメーションを作図
add_size <- 2
anim <- ggplot() + 
  geom_path(data = coord_circle_df, 
            mapping = aes(x = x, y = y, group = r), 
            color = "white") + # ノルム軸線
  geom_segment(data = coord_oblique_df, 
               mapping = aes(x = 0, y = 0, xend = x, yend = y, group = i), 
               color = "white") + # 角度軸線
  geom_segment(mapping = aes(x = c(-Inf, 0), y = c(0, -Inf), 
                             xend = c(Inf, 0), yend = c(0, Inf)), 
               arrow = arrow(length = unit(10, units = "pt"), ends = "last")) + # x・y軸線
  geom_text(data = coord_oblique_df, 
            mapping = aes(x = x, y = y, label = t_label, hjust = h, vjust = v), 
            parse = TRUE) + # 角度目盛ラベル
  geom_path(data = spiral_df, 
            mapping = aes(x = x, y = y, color = sgn_t_flag), 
            linewidth = 0.5) + # 螺旋
  geom_segment(data = angle_oblique_df,
               mapping = aes(x = 0, y = 0, xend = x, yend = y, linetype = type),
               linewidth = 1) + # 角度用の補助線
  geom_point(data = angle_point_df, 
             mapping = aes(x = x, y = y), 
             size = 5) + # 螺旋上の点
  geom_point(data = point_df, 
             mapping = aes(x = x, y = y, color = sgn_t_flag), 
             size = 4, shape = "circle open") + # 角度用の補助線上の点
  geom_text(data = angle_point_df, 
            mapping = aes(x = -Inf, y = Inf, label = var_label), 
            parse = TRUE, hjust = 0, vjust = -0.5) + # 変数ラベル
  gganimate::transition_manual(frames = frame_i) + # フレーム切替
  scale_color_manual(breaks = c("TRUE", "FALSE"), 
                     values = c("blue", "red"), 
                     labels = c(expression(theta >= 0), expression(theta < 0)), 
                     name = expression(sgn~theta)) + # 角度の正負
  scale_linetype_manual(breaks = c("origin", "negation"), 
                        values = c("solid", "dashed"), guide = "none") + # 角度の正負
  coord_fixed(ratio = 1, clip = "off", 
              xlim = c(-axis_size-add_size, axis_size+add_size), 
              ylim = c(-axis_size-add_size, axis_size+add_size)) + 
  labs(title = "Archimedes' spiral", 
       subtitle = "", # (変数ラベルの表示用)
       x = expression(x == r ~ cos~theta), 
       y = expression(y == r ~ sin~theta))

# gif画像を作成
gganimate::animate(
  plot = anim, 
  nframes = frame_num, fps = 10, 
  width = 600, height = 600, 
  renderer = gganimate::gifski_renderer()
)

 gganimate パッケージを利用してアニメーション(gif画像)を作成します。
 transition_manual() のフレーム制御の引数 frames にフレーム番号列 frame_i を指定します。
 animate()plot 引数にグラフオブジェクト、nframes 引数にフレーム数 frame_num を指定して、gif画像を作成します。また、fps 引数に1秒当たりのフレーム数を指定できます。

 フレームごとの \(\theta\) の点を黒色の点で示します。原点と \(\theta\) の点を通る半直線を黒色の線分、\(r\) が正の値と負の値を含む場合は反転した点を通る半直線を破線の線分で示します。
 x軸線の正の部分と半直線の偏角が変数 \(\theta\) に対応します。
 \(\theta + 2 i \pi\) の点( \(i\) を整数として \(\theta\) の前後 \(2 \pi\) 間隔の点)が線分上に並ぶのを確認できます。

パラメータと形状の関係

 続いて、定数(パラメータ)の値の変化に応じたグラフを作成して、パラメータと螺旋の形状の関係を確認します。

 フレーム数とパラメータを指定して、螺旋の描画用のデータフレームを作成します。

# フレーム数を指定
frame_num <- 101

# パラメータを指定
a_vals <- seq(from = -2, to = 2, length.out = frame_num)

# 周回数を指定
lap_num <- 8

# 螺旋の座標を作成
anim_spiral_df <- tidyr::expand_grid(
  frame_i = 1:frame_num, 
  t       = seq(from = -lap_num*pi, to = lap_num*pi, length.out = 1000)
) |> # フレームごとにラジアンを複製
  dplyr::mutate(
    a = a_vals[frame_i], 
    r = a * t, 
    x = r * cos(t), 
    y = r * sin(t)
  )
anim_spiral_df
## # A tibble: 101,000 × 6
##    frame_i     t     a     r     x        y
##      <int> <dbl> <dbl> <dbl> <dbl>    <dbl>
##  1       1 -25.1    -2  50.3  50.3 4.92e-14
##  2       1 -25.1    -2  50.2  50.1 2.52e+ 0
##  3       1 -25.0    -2  50.1  49.8 5.03e+ 0
##  4       1 -25.0    -2  50.0  49.4 7.51e+ 0
##  5       1 -24.9    -2  49.9  48.9 9.97e+ 0
##  6       1 -24.9    -2  49.8  48.2 1.24e+ 1
##  7       1 -24.8    -2  49.7  47.4 1.48e+ 1
##  8       1 -24.8    -2  49.6  46.5 1.71e+ 1
##  9       1 -24.7    -2  49.5  45.5 1.94e+ 1
## 10       1 -24.7    -2  49.4  44.4 2.16e+ 1
## # ℹ 100,990 more rows

 フレーム数 frame_num を指定して、frame_num 個の値を指定します。
 フレーム番号( 1 から frame_num までの整数)とラジアンの値の組み合わせを expand_grid() で作成します。パラメータの値ごとに螺旋用のラジアンを複製できます。

 目安として用いる螺旋上の点の描画用のデータフレームを作成します。

# 点数を指定
point_num <- lap_num * 4 + 1

# 螺旋上の点の座標を作成
anim_point_df <- tidyr::expand_grid(
  frame_i = 1:frame_num, 
  t       = seq(
    from = min(anim_spiral_df[["t"]]), 
    to   = max(anim_spiral_df[["t"]]), 
    length.out = point_num
  )
) |> # フレームごとにラジアンを複製
  dplyr::mutate(
    a = a_vals[frame_i], 
    r = a * t, 
    x = r * cos(t), 
    y = r * sin(t)
  )
anim_point_df
## # A tibble: 3,333 × 6
##    frame_i     t     a     r         x         y
##      <int> <dbl> <dbl> <dbl>     <dbl>     <dbl>
##  1       1 -25.1    -2  50.3  5.03e+ 1  4.92e-14
##  2       1 -23.6    -2  47.1 -1.27e-13  4.71e+ 1
##  3       1 -22.0    -2  44.0 -4.40e+ 1 -3.77e-14
##  4       1 -20.4    -2  40.8 -4.00e-14 -4.08e+ 1
##  5       1 -18.8    -2  37.7  3.77e+ 1  2.77e-14
##  6       1 -17.3    -2  34.6 -8.47e-14  3.46e+ 1
##  7       1 -15.7    -2  31.4 -3.14e+ 1 -1.92e-14
##  8       1 -14.1    -2  28.3  1.56e-14 -2.83e+ 1
##  9       1 -12.6    -2  25.1  2.51e+ 1  1.23e-14
## 10       1 -11.0    -2  22.0 -9.43e-15  2.20e+ 1
## # ℹ 3,323 more rows

 螺旋の対応関係の確認用に、変数(ラジアン)に関して一定間隔に点を描画することにします。
 点の数 point_num を指定し point_num 個のラジアンの値を作成して、螺旋の座標計算と同様に処理します。

 パラメータラベルの描画用のデータフレームを作成します。

# パラメータラベル用の文字列を作成
anim_label_df <- tibble::tibble(
  frame_i = 1:frame_num, # フレーム番号
  a       = a_vals, 
  param_label = paste0(
    "list(", 
    "r == a * theta", ", ", 
    "a == ", round(a, digits = 2), 
    ")"
  ) # パラメータラベル
)
anim_label_df
## # A tibble: 101 × 3
##    frame_i     a param_label                     
##      <int> <dbl> <chr>                           
##  1       1 -2    list(r == a * theta, a == -2)   
##  2       2 -1.96 list(r == a * theta, a == -1.96)
##  3       3 -1.92 list(r == a * theta, a == -1.92)
##  4       4 -1.88 list(r == a * theta, a == -1.88)
##  5       5 -1.84 list(r == a * theta, a == -1.84)
##  6       6 -1.8  list(r == a * theta, a == -1.8) 
##  7       7 -1.76 list(r == a * theta, a == -1.76)
##  8       8 -1.72 list(r == a * theta, a == -1.72)
##  9       9 -1.68 list(r == a * theta, a == -1.68)
## 10      10 -1.64 list(r == a * theta, a == -1.64)
## # ℹ 91 more rows

 フレームごとにパラメータラベル用の文字列を作成します。

 グラフサイズや軸線用の値を設定します。

# グラフサイズを設定
axis_size <- c(anim_spiral_df[["x"]], anim_spiral_df[["y"]]) |> 
  abs() |> 
  max() |> 
  ceiling()

# 目盛間隔を設定
step_val <- 12.5

# 軸線数を設定
circle_num <- axis_size %/% step_val
axis_size; circle_num
## [1] 51
## [1] 4

 グラフサイズ(の半分の値) axis_size、ノルム軸線の目盛間隔 step_val と軸線数 circle_num を設定します。この例では、ノルム軸の最大値によってグラフサイズを調整しません。

 「螺旋の作図」のときのコードで、ノルム軸線と角度軸線の描画用のデータフレームを作成します。

 螺旋のアニメーションを作成します。

# 螺旋のアニメーションを作図
anim <- ggplot() + 
  geom_path(data = coord_circle_df, 
            mapping = aes(x = x, y = y, group = r), 
            color = "white") + # ノルム軸線
  geom_segment(data = coord_oblique_df, 
               mapping = aes(x = 0, y = 0, xend = x, yend = y, group = i), 
               color = "white") + # 角度軸線
  geom_path(data = anim_spiral_df, 
            mapping = aes(x = x, y = y, color = t/pi), 
            linewidth = 1) + # 螺旋
  geom_point(data = anim_point_df, 
             mapping = aes(x = x, y = y, color = t/pi), 
             size = 4) + # 螺旋上の点
  geom_text(data = anim_label_df, 
            mapping = aes(x = -Inf, y = Inf, label = param_label), 
            parse = TRUE, hjust = 0, vjust = -0.5) + # パラメータラベル
  gganimate::transition_manual(frames = frame_i) + # フレーム切替
  coord_fixed(ratio = 1, clip = "off", 
              xlim = c(-axis_size, axis_size), 
              ylim = c(-axis_size, axis_size)) + 
  labs(title = "Archimedes' spiral", 
       subtitle = "", # (パラメータラベルの表示用)
       color = expression(frac(theta, pi)), 
       x = expression(x == r ~ cos~theta), 
       y = expression(y == r ~ sin~theta))

# gif画像を作成
gganimate::animate(
  plot = anim, 
  nframes = frame_num, fps = 10, 
  width = 600, height = 600, 
  renderer = gganimate::gifski_renderer()
)

 \(a\) の絶対値が大きいほど、\(r\) の絶対値が大きくなるため、螺旋のサイズが大きくなります。\(a\) の符号によってx軸・y軸に対して反転します。

アルキメデスの螺旋とsin曲線・cos曲線の関係

 最後は、アルキメデスの螺旋とsin関数・cos関数の曲線のグラフを作成して、変数(ラジアン)と座標や定数(パラメータ)と形状の関係を確認します。
 作図コードについては「archimedes_spiral.R」、作図の解説については「cos関数の可視化」などを参照してください。

変数と座標の関係

 変数に応じて移動するアルキメデスの螺旋・sin曲線・cos曲線上の点のアニメーションを作成します。

 螺旋上の点(黒色の点)のx軸方向の変化(青色の点)とy軸方向の変化(赤色の点)が、それぞれcos関数曲線上の点とsin関数曲線上の点に対応しているのを確認できます。
 \(\theta\) の絶対値が大きいほど、\(r\) の絶対値が大きくなるため、曲線の振幅が大きくなります。

パラメータと形状の関係

 パラメータに応じて変化するアルキメデスの螺旋・sin曲線・cos曲線のアニメーションを作成します。

 対応関係の目安として、\(\theta\) に応じて曲線を色付けし、また一定間隔で点を打っています。 螺旋のx軸方向の形状とy軸方向の形状が、それぞれcos関数曲線の形状とsin関数曲線の形状に対応しているのを確認できます。
 \(a\) の絶対値が大きいほど、\(r\) の絶対値が大きくなるため、曲線の振幅が大きくなります。\(a\) の符号によってそれぞれの軸に対して反転します。

 この記事では、アルキメデスの螺旋のグラフを作成しました。次の記事では、対数螺旋のグラフを作成します。